function rest=maxL2err(u,term,errtype)

n=floor((size(u,1)+term+0.25)/term);    h=1/n;   unimesh=(0:1/n:1)';

if(strcmp(errtype,'max'))
    utmp=[u0ex(unimesh(2:n));u1ex(unimesh(2:n));u2ex(unimesh(2:n));u3ex(unimesh(2:n))];
    uex=utmp(1:term*(n-1));
    rest=max(abs(uex-u));
elseif(strcmp(errtype,'L2'))
    load hermitdata;
    
switch term
    case 1
        u0p=u(1:n-1);         
    case 2
        u0p=u(1:n-1); u1p=u(n:2*n-2);
    case 3
        u0p=u(1:n-1); u1p=u(n:2*n-2);   u2p=u(2*n-1:3*n-3); 
        
    case 4
        u0p=u(1:n-1); u1p=u(n:2*n-2);   u2p=u(2*n-1:3*n-3);   u3p=u(3*n-2:4*n-4);
        err = zeros(n-1,1);
        for ii=1:n-1
            x=ii*h;
            errorhd1=@(z) (log(1+z*x)./log(1+z) ...
                 - u0p(ii) - u1p(ii)*z/sig - u2p(ii)/sqrt(2)*(z.^2/sig^2-1) ...
                 - u3p(ii)/sqrt(6)*(z.^3/sig^3-3*z/sig)).^2;
            errorhd2=@(z) ((1+z)*x/(1+z*x) ...
                 - u0p(ii) - u1p(ii)*z/sig - u2p(ii)/sqrt(2)*(z.^2/sig^2-1) ...
                 - u3p(ii)/sqrt(6)*(z.^3/sig^3-3*z/sig)).^2;
            tmphd1=@(z) errorhd1(z).*gaussian(z);
            tmphd2=@(z) errorhd2(z).*gaussian(z);
            epsilon=10^-4;
            err(ii)=quadl(tmphd1,-1,-epsilon,10^-10);%+quadl(tmphd2,-epsilon,epsilon,10^-10)+quadl(tmphd1,epsilon,1,10^-10);
        end
        rest=sqrt(max(err));        
end
end